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We study phase entrainment of Kuramoto oscillators under different conditions on the interaction 
range and the natural frequencies. In the first part the oscillators are entrained by a pacemaker 
acting like an impurity or a defect. We analytically derive the entrainment frequency for arbitrary 
interaction range and the entrainment threshold for all-to-all couplings. For intermediate couplings 
our numerical results show a reentrance of the synchronization transition as a function of the coupling 
range. The origin of this reentrance can be traced back to the normalization of the coupling strength. 
In the second part we consider a system of oscillators with an initial gradient in their natural 
frequencies, extended over a one-dimensional chain or a two-dimensional lattice. Here it is the 
oscillator with the highest natural frequency that becomes the pacemaker of the ensemble, sending 
out circular waves in oscillator-phase space. No asymmetric coupling between the oscillators is 
needed for this dynamical induction of the pacemaker property nor need it be distinguished by a 
gap in the natural frequency. 
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I. INTRODUCTION 

Synchronization in the sense of coordinated behavior in 
time is essential for any efficient organization of systems, 
natural as well as artificial ones. In artificial systems 
like factories the sequence of production processes should 
be synchronized in a way that time-and space consum- 
ing storage of input or intermediate products is avoided. 
The same applies to natural systems like networks of cells 
which obviously perform very efficiently in fulfilling a va- 
riety of functions and tasks. In a more special sense, the 
synchronized behavior refers to oscillators with almost 
identical individual units, in particular to phase oscil- 
lators with continuous interactions as described by the 
Kuramoto model These sets of limit-cycle oscilla- 
tors describe synchronization phenomena in a wide range 
of applications Q. One of these applications is pattern 
formation in chemical oscillatory systems 0, described 
by reaction-diffusion systems. A well-known example for 
such a chemical system is the Belousov-Zabotinsky re- 
action, a mixture of bromate, bromomalonic acid and 
ferroin, which periodically changes its color correspond- 
ing to oscillating concentrations of the red, reduced state 
and the blue, oxidized state. In these systems, expand- 
ing target-like waves and rotating spiral waves are typical 
patterns 0. As Kuramoto showed (as e.g. in pj), these 
dynamical systems can be well approximated by phase os- 
cillators if the interaction is weak. For reaction-diffusion 
equations with a nonlinear interaction term he predicted 
circular waves with strong similarities to the experimen- 
tal observation of target patterns. In these systems, de- 
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fects or impurities seem to play the role of pacemakers, 
driving the system into a synchronized state. Therefore 
he treated pacemakers as local "defect" terms leading to 
heterogeneities in the reaction-diffusion equations. 

In this paper we consider two types of pacemakers. 
In the first part they are introduced as "defects" in the 
sense that they have a different natural frequency from 
the rest of the system, whose oscillators have either ex- 
actly the same frequency or small random fluctuations 
about some common average by assumption, where in 
the second case the pacemaker's frequency is clearly dif- 
ferent from the values of a typical fluctuation. Actually, 
we will not consider the latter case, since small fluctua- 
tions do not lead to any qualitative change. The ad-hoc 
distinction of the pacemaker may reflect natural and ar- 
tificial systems with built-in impurities like those in the 
Belousov-Zhabotinsky system. We have recently stud- 
ied these systems for nearest neighbor interactions 
and regular networks in d dimensions for periodic and 
open boundary conditions. For such systems we have 
shown that it is the mean distance of all other nodes 
from the pacemaker (the so-called depth of the network) 
that determines its synchronizability. Periodic bound- 
ary conditions facilitate synchronization as compared to 
open boundary conditions. Here we extend these results 
to larger interaction range. The effect of an extended 
interaction range on oscillatory systems without pace- 
maker were considered in 6] . We analytically derive the 
entrainment frequency for arbitrary interaction range be- 
tween next-neighbor and all-to all interactions. We de- 
rive the entrainment window for all- to-all couplings an- 
alytically, for intermediate couplings numerically. The 
entrainment window depends non-monotonically on the 
interaction range, so that the synchronization transition 
is reentrant. In our system the reentrance can be ex- 
plained in terms of the weight of the interaction term 
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and the normalization of the coupling strength. 

In the second part we consider a system of oscillators 
without a " defect" , all oscillators on the same footing up 
to the difference in their natural frequencies. An obvious 
choice for the natural frequencies would be a Gaussian 
distribution or another random distribution to describe 
fluctuations in natural frequencies in otherwise homoge- 
neous systems, without impurities. Such systems have 
been studied by Blasius and Tonjes for a Gaussian 
distribution. The authors showed that in this case an 
asymmetric interaction term is needed to have synchro- 
nization, this time driven by a dynamically established 
pacemaker, acting as a source of concentric waves. This 
result shows that conditions exist under which a system 
establishes its own pacemaker in a "self-organized" way. 
In contrast to the random distribution of natural frequen- 
cies we are interested in a (deterministic) gradient distri- 
bution, in which the natural frequencies linearly decrease 
over a certain region in space. Such a deterministic gra- 
dient alone is probably not realistic for natural systems, 
but considered as a subsystem of a larger set of oscilla- 
tors with natural frequency fluctuations it may be real- 
ized over a certain region in space. Here we suppress the 
fluctuations and focus on the effect of the gradient alone. 
As it turns out, the asymmetry in the natural frequency 
distribution is sufficient for creating a pacemaker as cen- 
ter of circular waves in oscillator-phase space, without 
the need for an asymmetric term in the interaction. The 
oscillator with the highest natural frequency becomes the 
source of circular waves in the synchronized system. We 
call it dynamically induced as its only inherent difference 
is its local maximum at the boundary in an otherwise 
"smooth" frequency distribution. The type of pattern 
created by the pacemaker has the familiar form of cir- 
cular waves. Beyond a critical slope in the natural fre- 
quencies, full synchronization is lost. It is first replaced 
by partial synchronization patterns with bifurcation in 
the frequencies of synchronized clusters, before it gets 
completely lost for too steep slopes. 

The outline of the paper is as follows. In the first 
part we treat the pacemaker as defect and derive the 
common entrainment frequency for arbitrary interaction 
range and topology (section [n|. Next we determine the 
entrainment window, analytically for all-to all coupling 
and numerically for intermediate interaction range (sec- 
tioning Here we see the reentrance of the transition as a 
function of the interaction range (section [HI). In the sec- 
ond part (section IIII|) we consider dynamically induced 
pacemakers, first without asymmetric interaction term, 
for which we analytically derive the synchronization tran- 
sition as a function of the gradient in the frequencies 
C sections IIII All . In section All Bl we add an asymmetric 
term, so that the pattern formation is no longer surpris- 
ing due to the results of [?J , but the intermediate patterns 
of partial synchronization are different. In section HVl we 
summarize our results and conclusions. 



II. PACEMAKER AS DEFECT IN THE SYSTEM 

The system is defined on a network, regular or small- 
world like. To each node i, i — 0, . . . ,N, we assign a 
limit-cycle oscillator, characterized merely by its phase 
ipi, which follows the dynamics 

is 

(pi = cj + 6 tiS Auj + — y] Ajj sin (ipj - ipi) , (1) 

ki . 

with the following notations. The frequency lu denotes 
the natural frequencies of the system. In this system we 
treat the pacemaker as a defect. It is labelled by s and has 
a natural frequency that differs by Aw from the frequency 
of the other oscillators having all the same frequency v. 
8i j denotes the Kronecker delta. The constant K > 
parameterizes the coupling strength. We consider regular 
networks and choose 

A 0,i — r j,i ' 

Tj t i the distance between nodes j and i, that is 

r Sti = wm([\j-i\ , (N + l)-\j-i\]) (2) 

on a one-dimensional lattice with periodic boundary con- 
ditions. In two or higher dimensions it is the shortest 
distance in lattice links. The parameter < a < oo 
tunes the interaction range. Alternatively we consider 
Aji as the adjacency matrix on a small-world topology: 
Aj i = I if the nodes j and i are connected and Aj i — 
otherwise). Moreover, ki = Y^j A?>* is the degree of the 
i-th node, it gives the total number of connections of this 
node in the network. This system was considered before 
in |5| for nearest-neighbor interactions (a — > oo) on a 
d-dimensional hypercubic lattice and on a Cayley tree. 
Here we extend the results to long-range interactions via 
< a < oo. 



Entrainment frequency and entrainment window 

In the appendix of [j| we derived the common entrain- 
ment frequency Q in the phase-locked state, for which 
(pi = for all i, to be given as 

n = a w A- (3) 

in the rotated frame, in which the natural frequency u) is 
zero. In particular such a result is obtained directly from 
system applying the transformation to the rotating 
frame ifi — > ifi + Lot Vz, multiplying the z-th equation by 
hi, summing over all equations and then using the fact 
that Yli^2j Aj,i sin ((fj — ipi) = 0, because the symme- 
try of the adjacency matrix A and the antisymmetry of 
the sin-function. Whenever fej = k is independent of i 
as it happens for periodic boundary conditions and ho- 
mogeneous degree distributions, in particular for a = 
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(all-to-all coupling), the common frequency (in the non- 
rotated frame) is given by 



n 



Auj 
N+1 



independent on a. In the limit TV — » oo the system syn- 
chronizes at lo. 

Analytical results for the entrainment window were de- 
rived for open and periodic boundary conditions and 
a — ► oo in For a = 0, the result is the same as 
for a — > oo, periodic boundary conditions and d > 1 



Auj 



1 

TV+ 1 



(4) 



This is seen as follows. Since for all-to-all couplings the 
mean-field approximation becomes exact, it is natural to 
rewrite the dynamical equations Q in terms of the order 
parameter Re % ^ defined according to 

ft e # = —X^ e lip > . 
TV ^ 

This quantity differs from the usual order parameter just 
by the fact that the sum runs over j ^ s. Eq.sQ then 
take the form 

fl = Auj + K R sm(ijj — (p s ) , for i = s 

£1 = |$ sin(</? s — (fit ) + K R sxa.{ijj — ipi) , V i ^= s ' 

(5) 

Using ifi = i/j for all i 7^ s as one of the possible so- 
lutions and the fact that the sin-function is bounded, 
these equations imply the entrainment window as given 
by (@J with TV the total number of nodes minus 1. As 
we have shown in the same dependence in terms of 
TV holds for a d-dimensional hypercubic lattice with TV' 
nodes per side and nearest neighbor couplings, for which 
TV = (TV' + l) d — 1, so that the entrainment window for 
nearest neighbor interactions and on a regular lattice de- 
cays exponentially with the depth of the network in the 
limit of infinite dimensions d ; for random and small- 
world networks the same type of decay was derived in 

a 

The intermediate interaction range < a < 00 can- 
not be treated analytically. Here the numerical inte- 
gration of the equations (IJ (via a fourth-order Rungc- 
Kutta algorithm with time step dt = 0.05 and homo- 
geneous initial conditions) shows that the entrainment 
window is even smaller than in the limiting cases a = 
or a — 00 for otherwise unchanged parameters (that is 
fixed size TV, dimension d, coupling K). For these cases 
formula @ gives an upper bound on the entrainment 
window, as one can easily see from Eq^ for i = s be- 
cause T,j^s sin {fj ~ ¥s)/rls < Ei^V'ts = k s > Va. 
Actually the entrainment window becomes smaller for 
intermediate interaction range, a feature leading to the 
reentrance of the synchronization transition as function 
of a, as we shall see in the figures ^ and [3 below. 



Reentrance of the phase transition 

Reentrance phenomena are observed in a variety of 
phase transitions, ranging from superconductor- insulator 
transitions as function of temperature and noise- 
induced transitions to chaotic coupled-map lattices 
as function of the coupling ,11.]. In |T^ | a reentrance 
phenomenon is discussed as an artifact of the approx- 
imation. Reentrance of phase transitions is challeng- 
ing as long as it appears counterintuitive. For exam- 
ple, it is counterintuitive when synchronization is first 
facilitated for increasing coupling and later gets blocked 
when the coupling exceeds a certain threshold, so that 
synchronization depends non-monotonically on the cou- 
pling. As we see in Fig. nf° r a one-dimensional system 
with periodic boundary conditions and for various sys- 
tem sizes TV, the entrainment window \Au/K\ c depends 
non-monotonically on the interaction range, parameter- 
ized by a. This dependence is easily explained by looking 
at the total weight coming from the sum of interactions. 
This contribution depends non-monotonically on a, as 
E • sin(ipj — 4>i)/ r f,i decreases for increasing a, while the 
coupling strength increases under the same change of a, 
so that there is a competition between the two factors 
leading to a minimum of the weight for a certain value 
of a. In contrast, if we scale the coupling strength by 
a constant factor K/N, independently of the degree of 
the node i, the overall weight of the sum decreases for 
increasing a and along with this, the drive to the syn- 
chronized state. As it is shown in the inset of Fig. 2] the 
threshold monotonically decreases with a if the coupling- 
strength is normalized by 1/TV independently of the ac- 
tual degree, as expected from our arguments above. The 
second inset shows the TV-dependence of the value of a 
for which synchronization is most difficult to achieve, it 
converges to a — 3 for large TV. Similar non-monotonic 
behavior is observed if the interaction range is tuned via 
the number k of nearest neighbors on a ring rather than 
via a, see Fig. |2K an d for a small- world topology [l^, see 
Fig. [33. Here the small- world is constructed by adding 
random shortcuts with probability p from a regular ring 
with nearest neighbor interactions (p = 0) to an all-to-all 
topology (p= 1) 0. 



III. DYNAMICALLY INDUCED PACEMAKERS 

In view of dynamically induced pacemakers we study 
the dynamical system 



tfi = t*)j + — ^ Aj,iT((fij - ipi) 



(0) 



without a pacemaker. The interaction term T is given by 

T(ipj - ifi) = sin(ipj-ipi) + y[l - cos (cpj - ipi)] , 

(7) 

where 7 = corresponds to antisymmetric interaction, 
7 > to asymmetric interaction. In the following we 
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Figure 1: Normalized entrainment window as function of the 
interaction range parameterized via a, for a one- dimensional 
lattice with N + 1 sites, periodic boundary conditions and for 
interaction strength between the oscillators i and j suppressed 
by the power a of their distance r^i. The main plot indicates 
the reentrance of the synchronization transition. The upper 
inset shows the size dependence of a m , for which the entrain- 
ment window becomes minimal, it saturates at a ~ 3. The 
lower inset shows the monotonic dependence on a for constant 
coupling strength K/N . The figure was obtained for homoge- 
neous initial conditions. The qualitative features remain the 
same for random initial conditions. 



consider the cases 7 = and 7 = 2. The value 7 = 2 is 
chosen for convenience in order of comparing our results 
with those of Blasius and Tonjes Q , but different values 
of 7 > lead to qualitatively the same features. 
We focus on lattices with periodic boundary conditions 
in one and two dimensions. Moreover we consider only 
nearest neighbor interactions, so that in d = 1 we have 
ki = 2 Vi and fej = 4 Vi in d = 2. The natural frequencies 
u>i are chosen from a gradient distribution 



Ujj 



Alu 



Alu 



max, n 



(8) 



where rij is the network distance between the oscilla- 
tors i and j, given by J2J in one dimension and by the 
shortest path along the edges of the square lattice in two 
dimensions. The oscillator s is the oscillator with the 
highest natural frequency. For a one-dimensional lattice 
with periodic boundary conditions (ring) and N + 1 os- 
cillators, choosing s = for simplicity, max^ 770 = N r , 
with N' = N/2 for N even and N' = (N + l)/2 for N 
odd. In the following we choose N as even; the case of 
N odd is very similar. 



A. Pattern formation for 7 = 

1. Entrainment frequency 
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Figure 2: Normalized entrainment window as function of the 
number of neighbors k (A) and of the probability p for adding 
shortcuts in a small-world topology (B) for a one-dimensional 
lattice with periodic boundary conditions and N + 1 sites. 
Average values in (B) are taken over 50 realizations, error 
bars correspond to the standard deviation. Again the non- 
monotonic behavior indicates reentrance of the synchroniza- 
tion transition. 



using the antisymmetry of T for 7 = 0, we obtain 



(N + l)n = = Aw + 2Aw Ei=i Z (1 - 2*/JV) 



N/2 , 



Alu 



i + n-a/n mamm 



= Acf 



leading to 



fl = Alj- 



N 



2(N + 1) ' 
so that — > Aw/2 for N — * 00. 

2. Entrainment window 



(9) 



Next we derive the common frequency in case of phase- 
locking tfii = fi , \fi. Summing the equations of @ and 



It is convenient to rewrite the differential equations © 
in terms of phase lags Qi = Lpi — ipi—x, so that the equation 
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for the oscillator at position reads 
ft = u) + if sin (0i ) sin (0i ) 
For the other oscillators we obtain 



K 



Q = U! + K/2 [sin (0 2 ) - sin (0 X )] 
=> sin(0 2 ) = 2^^+ sin (0 X ) 



n = lu 2 + K/2 [sin (0 3 ) - sin (0 2 )] 
=>■ sin (0 3 ) 



2^+sin(0 2 ) 



so that in general 

Q = Ui-i + K/2 [sin (0,) - sin (0;_i)] 
=> sin(fli) = 2^^+sin(0 4 _ 1 ) 



K 



The former sum can be performed as 



E l—l n ! 
.7=1 Z " 



A" 



A" 



E!C 1 i[«-Ao;(l-2i/7V)] 



£ (n - a^) (*-i) + Aw-|i^ 



AT 



from which 

sin(0i) 
or, using Eq.®, 



fi-Aw , , Aw 2i 2 - 2i 
(2i-l) + 



K 



sin 



= ^(2i-l) 



AT 



A 



(1 - 2i) 



7V+ 2 
2(A+1) 



2 (A+1 
. 2i(i-l) 



- 1 



N 



Auj 2i(i — 1) 
A A 



A 



(10) 

For given N, we first determine the position i m at which 
the right-hand-side of Eq. (|10|l takes its minimum (maxi- 
mum) value as a function of i. This value is constrained 
by the left-hand-side of (fTT)f) to be out of [—1,1]. For i 
treated as a continuous variable the derivative of Eq. i|l(J|) 
with respect to i yields 



N(N + 2) 1 



4(iV+l) 2 ' 

On the other hand 

i m > 1 & N > V2 , 

and 

N r- 

im<— iV>l + V3, 



(11) 



(12) 



(13) 



so that for N 1 , 1 < i m < N/2 is certainly satisfied. 
Upon inserting i m from (|ll|l we have 
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Figure 3: Critical threshold for a one-dimensional lattice as 
function of its linear size N, full line corresponds to the ana- 
lytic prediction, full dots to numerical results for N > 2, the 
symmetry breaking parameter is chosen as 7 = 0. 



The critical threshold is then given by 



Aw 
~K 



8N(N+ l) 2 



N 2 (N + 2) 2 + A(N + 1) 2 



(14) 



This relation is plotted in Fig. |3 where the full line rep- 
resents the theoretical prediction, while the full dots cor- 
respond to the threshold values, numerically determined 
by integrating the system of differential Eq.s ©, starting 
with N = 2. (The case of N = 2 violates the constraint 
(|13|l in agreement with the fact that for N = 2 the sys- 
tem degenerates to the former case of having a pacemaker 
as a defect and two oscillators with natural frequencies 
u>x = w 2 = 0.) Furthermore it should be noticed that 
Eci. (|14|) scales as 1/N for large values of N, the large- 
TV behavior is, however, not visible for the range of N, 
(N < 20), plotted in Fig. 

For a ring with maximal distance N' = N/2 and given 
N and K, the slope Aw/iV' in the natural frequencies de- 
termines whether synchronization is possible or not. The 
oscillators synchronize if the values of iV lead to a ratio 
I Auj/K I below the bound \ Au/K\n of Ea.l(T4'jl. named as 
in section [n] Eq.@, but with Auj now standing for the 
maximal difference between the highest and the lowest 
natural frequencies. If for given N the slope, parame- 
terized by Aw, exceeds this threshold, no entrainment is 
possible. For N — > 00 the allowed slope goes to zero, the 
system does no longer synchronize. To further character- 
ize the synchronization patterns, we measure the station- 
ary frequency fi of the individual oscillators i, defined 
according to 



fi = Ihn 



<Pi{t + t ) - ipi(t Q ) 
t-t 



(15) 



well as its average / = jA^ 2~2i=o fi an d variance 

1/2 

. The average frequency 



N+ 



+ 1 2~2i=0 (fi f ) 
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and variance are introduced to characterize the interme- 
diate patterns of partial synchronization and to distin- 
guish it from the case 7 = 2 in section IIII Bl In all 
numerical calculations we set from now on K — 1. Fig. 
0] shows the natural frequencies (u>i) (represented as full 
line) and the stationary (/,) frequencies (shown as dots) 
for four slopes Aw = 0.01,0.1,0.2,0.3, respectively, as 
function of the position on a linear chain of 100 oscilla- 
tors, of which we show only one half due to the symmet- 
ric arrangement. The case of Aw = 0.01 corresponds to 
full synchronization, while the other three figures repre- 
sent partial synchronization of two (Aw = 0.1) and more 
clusters. It should be noticed that the two clusters for 
(Aw = 0.1) would violate the condition iflljl if they were 
isolated clusters of 50 oscillators each with a slope deter- 
mined by Aw = 0.1, but due to the nonlinear coupling 
to the oscillators of the second cluster, the partial syn- 
chronization is a stationary pattern. Fig. displays the 
frequencies as function of the slope, again parameter- 
ized by Auj. Here we have chosen N = 20. We clearly see 
the bifurcation in frequency space starting at the critical 
value of Alu ~ 0.36 and ending with complete desyn- 
chronization, in which the stationary frequencies equal 
the natural frequencies. For frequency synchronization 
the phases are locked, but their differences to the phase 
of say the oscillator s = (ipi — tfi — tpo) increases non- 
linear ly with their distance from s, as seen in Fig. [S] 
Therefore the distance between points of the same ip de- 
creases with i. 

In two dimensions we simulated 100 x 100 oscillators on 
a square lattice with periodic boundary conditions and 
the oscillator s with the highest natural frequency placed 
at the center of the square lattice. The natural frequen- 
cies of the oscillators are still given by Eq.JSJ). Their 
spatial distribution looks like a square pyramid centered 
at the middle of the square lattice. Fig. [TJdisplays sin?/*, 
on this grid. It exhibits stationary patterns after 2 • 10 4 
steps of integration, for six choices of Auj, all above the 
synchronization threshold. We see some remnants of syn- 
chronization, most pronounced in Fig. in which a cir- 
cular wave is created at the center at s, and coexists with 
waves absorbed by sinks at the four corners of the lattice. 
The projection of Fig. [7JA on one dimension corresponds 
to Fig. with a bifurcation into two cluster-frequencies. 
Note that it is again the oscillator with the highest nat- 
ural frequency that becomes the center of the outgoing 
wave, while the corners with the lowest natural frequen- 
cies (w; — 0) become sinks. Although Fig. [7p shows 
almost no remnant of synchronization, it is interesting to 
follow the time evolution towards this "disordered" pat- 
tern via a number of snapshots, as displayed in Fig. [HI af- 
ter 9 [22.5] (A), 41 [102.5] (B) , 150 [375] (C) , 393 [982.5] 
(D), 490 [1225] (E), 1995 [4987.5] (F) integration steps 
T [T Auj time steps in units of natural frequency(Aw) 
of Eq|S]. While the pattern of (A) would be stationary 
in case of full synchronization, here it evolves after it- 
erated reflections to that of (F). The pattern of (F) is 
stationary in the sense that it stays "disordered", on a 
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Figure 4: Natural frequencies uji (full lines) and measured, 
stationary frequencies fi (full dots) as function of the position 
i of the oscillators on a one-dimensional lattice, for N = 100, 
periodic boundary conditions, and four values of the maximal 
natural frequency Alu. The oscillator with maximal natural 
frequency Aw is at position s = 0. Only half of the oscil- 
lators up to index 50 are shown because of the symmetric 
arrangement. The symmetry breaking parameter is chosen as 
7 = 0. 
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Figure 5: Bifurcation in cluster-frequency space above a crit- 
ical slope in the natural frequencies, parameterized by Aw, 
the maximal natural frequency. Numerically measured fre- 
quencies fi (normalized such that fi/ Auj £ [0, 1]) as function 
of Aw for a one-dimensional lattice with N = 20 and periodic 
boundary conditions. The symmetry breaking parameter is 
chosen as 7 = 0. 



"microscopic" scale it shows fluctuations in the phases. 
The evolution takes a number of time steps, since the in- 
teraction is mediated only via nearest neighbors and not 
of the mean-field type (all-to-all coupling) . 
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Figure 6: Phase differences ipi = ifi — ipo (from the oscilla- 
tor at i = 0) as function of the oscillator index i for a one- 
dimensional lattice, with TV = 100 and periodic boundary 
conditions, for different values of the maximal natural fre- 
quency Acj and 7 = 0. The difference depends nonlinear ly 
on the distance from the oscillator at i = 0. The "steps" in 
these plots are due to the projection of the variables ipi onto 
the interval [0, 2n). 



Figure 8: Time evolution of the "disorder" in Fig. [Tf 1 for 
symmetry breaking parameter 7 = 0. All plots are realized 
for Au; = 50, after a different number of integration steps 
[time in natural units of Aw]: A) 9 [22.5] B) 41 [102.5] , C) 
150 [375] , D) 393 [982.5] , E) 490 [1225], F) 1995 [4987.5]. 
Only the pattern of F) is stationary on a global scale, i.e., up 
to local fluctuations in the phases. 
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Figure 7: Pattern formation for symmetry breaking param- 
eter 7 = on a two-dimensional lattice with TV = 100 os- 
cillators per side. The oscillator s with the maximal natural 
frequency is placed at the center the lattice. We plot sinipi 
for every oscillator i. The colors vary between white, corre- 
sponding to —1, and black, corresponding to 1. The plots are 
realized after 2 ■ 10 4 integration steps and for different values 
of Aw [therefore after different time in natural units of Aw]: 
A) 0.1 [10 2 ] , B) 0.5 [5-10 2 ] , C) 1.0 [10 3 ] , D) 5.0 [5-10 3 ] , E) 
10 [10 4 ] , F) 50 [5 ■ 10 4 ]. A cross section through the pattern 
of A) would reveal a frequency distribution as in Fig|4*B. 



B. Pattern formation for 7 = 2 

Next let us consider system (JSJ) for 7 — 2, for which 
we expect pattern formation from the results of due 
to the broken antisymmetry of sin via a cos-term in the 
interaction. The main difference shows up in the par- 



tial synchronization patterns above the synchronization 
threshold. As it is evident from Fig.s 191 and 1101 the size 
of the one synchronized cluster, filling the whole lattice 
below the threshold, shrinks with increasing Aw, since 
the oscillators with the highest natural frequencies de- 
couple from the cluster, but the remaining oscillators do 
not organize in other synchronized clusters, the cluster- 
frequency does no longer bifurcate as before. Along with 
this, the two-dimensional stationary patterns change to 
those of Fig. El whereas the evolution towards Fig. I11F 
is now displayed in Fig.s^l w ith snapshots taken after 9 
[22.5] (A), 48 [120] (B), 160 [400] (C), 277 [692.5] (D), 386 
[965] (E), and 2926 [7315] (F) integration steps [T Aw as 
above] . 

A further manifestation of the difference in the partial 
synchronization patterns is seen in the average frequency 
/ and the variance a as function of Aw (Fig. I13fl . Above 
the synchronization transition the average frequency re- 
mains constant for 7 = 0, but increases for 7 = 2 in 
agreement with Fig.s and 1101 the variance increases 
in both cases (7 = and 7 = 2) above the transition, 
so that it may serve as order parameter. Therefore we 
can tune the synchronization features via the slope of the 
gradient. 

The main qualitative feature, however, is in common to 
both systems with and without antisymmetric coupling: 
the oscillator with the highest natural frequency becomes 
the center of outgoing circular waves, it is dynamically 
established as pacemaker. The patterns, seen here in 
case of full and partial synchronization, are quite similar 
to those predicted by Kuramoto [lj for reaction-diffusion 
systems with continuous diffusion terms, and to those 
experimentally observed in chemical oscillatory systems 
i- 




Figure 9: Natural and measured frequencies for symmetry 
breaking parameter 7 = 2. The frequency of the synchro- 
nized cluster decreases and the cluster size shrinks, until all 
oscillators keep their natural frequency in a completely desyn- 
chronized state. 




Figure 11: Pattern formation for symmetry breaking param- 
eter 7 = 2. We use the same parameter choice as in Fig. 
|7| Now Fig.s A), B), C) correspond to full synchronization, 
while one-dimensional cross sections through D), E), and F) 
would reveal distributions of /j as in Fig|5]C)-D). 
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Figure 10: Bifurcation in cluster-frequency for symmetry 
breaking parameter 7 = 2. Above the synchronization thresh- 
old the size of the original synchronized cluster shrinks with 
Aw. Oscillators outside this cluster stay isolated in contrast 

to Fig. oa 
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Figure 12: Time evolution for symmetry breaking parameter 
7 = 2, otherwise same as in Fig. |H| This time the plots 
are realized after 9 [22.5] (A), 48 [120] (B), 160 [400] (C), 
277 [692.5] (D), 386 [965] (E), and 2926 [7315] (F) time steps 
[time in units of Aw]. Only the disordered pattern of (F) is 
stationary up to local fluctuations in the phases. 



IV. SUMMARY AND CONCLUSIONS 

For pacemakers implemented as defect we have ex- 
tended former results on the entrainment frequency and 
the entrainment window to arbitrary interaction range, 
analytically for a = 0, numerically for intermediate 
< a < 00. For large dimensions the entrainment 
window decays exponentially with the average distance 
of nodes from the pacemaker so that only shallow net- 
works allow entrainment. The synchronization transition 
is reentrant as function of a (or k, the number of neigh- 
bors on a ring, or p, the probability to add a random 
shortcut to the ring topology). The entrainment gets 
most difficult for a ~ 3 and large N, while its most easily 



achieved for next-neighbor and all-to-all couplings. This 
reentrance is easily explained in terms of the normaliza- 
tion of the coupling strength. For the same system with- 
out a pacemaker, but with a gradient in the initial natural 
frequency distribution the oscillator with the highest nat- 
ural frequency becomes the center of circular waves, its 
role as a pacemaker is dynamically induced without the 
need for an asymmetric term in the interaction. Here we 
analytically determined the synchronization transition as 
function of the gradient slope. Above some threshold, full 
synchronization on one or two-dimensional lattices is lost 
and replaced by partial synchronization patterns. These 
patterns depend on the asymmetry parameter 7. For 
7 = we observe a bifurcation in frequency space, for 
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7 = 2 the one synchronized cluster shrinks in its size, be- 
fore for even steeper slopes synchronization is completely 
lost. For artificial networks these results may be used to 
optimize the placement and the number of pacemakers if 
full synchronization is needed, or to control synchroniza- 
tion by tuning the slope of natural frequency gradients. 



Figure 13: Average frequency and variance as defined in the 
text (both divided by Aoj) as function of the maximal nat- 
ural frequency Acu, in case of a one-dimensional lattice with 
N = 100 oscillators and for values of the symmetry breaking 
parameter 7 = and 7 = 2. Fhe raise in the variance from 
zero indicates the transition to the desynchronized phase. 
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